#include "EquationSolver.h"
#define pi acos(-1)
const double l = 89, h = 49, beta1 = 11.5/180*pi;

//题目F要求：In the design of all-terrain vehicles, it is necessary to consider the failure of the vehicle when attempting to negotiate two types of obstacles. One type of failure is called hang-up failure and occurs when the vehicle attempts to cross an obstacle that causes the bottom of the vehicle to touch the ground. The other type of failure is called nose-in failure and occurs when the vehicle descends into a ditch and its nose touches the ground.
//The above figure shows the components associated with the nose-in failure of a vehicle. The maximum angle α that can be negotiated by a vehicle when β is the maximum angle at which hang-up failure does not occur satisfies the equation A sin α cos α + B sin2 α − C cos α − E sin α = 0,where A = lsin β1, B = l cos β1, C = (h + 0.5D) sin β1 − 0.5D tan β1, E = (h + 0.5D) cos β1 − 0.5D.
//(a) Use Newton’s method to verify α ≈ 33◦ when l = 89 in., h = 49 in., D = 55 in. and β1 = 11.5◦.
//(b) Use Newton’s method to find α with the initial guess 33◦for the situation when l, h, β1 are the same as in part (a) but D = 30 in..
//(c) Use the secant method (with another initial value as far away as possbile from 33◦) to find α. Show that you get a different result if the initial value is too far away from 33◦; discuss the reasons.

class Fun1 : public Function
{
public:
    double operator()(double _x)
    {
        double y;
        double D = 55;
        double A = l*sin(beta1);
        double B = l*cos(beta1);
        double C = (h + 0.5*D)*sin(beta1) - 0.5*D*tan(beta1);
        double E = (h + 0.5*D)*cos(beta1) - 0.5*D;
        y = A*sin(_x)*cos(_x) + B*pow(sin(_x),2) - C*cos(_x) - E*sin(_x);
        return y;
    }
    double diff(double _x)
    {
        double y;
        double D = 55;
        double A = l*sin(beta1);
        double B = l*cos(beta1);
        double C = (h + 0.5*D)*sin(beta1) - 0.5*D*tan(beta1);
        double E = (h + 0.5*D)*cos(beta1) - 0.5*D;
        y = -A*pow(sin(_x),2) + A*pow(cos(_x),2) + 2*B*sin(_x)*cos(_x) + C*sin(_x) - E*cos(_x);
        return y;
    }
};

class Fun2 : public Function
{
public:
    double operator()(double _x)
    {
        double y;
        double D = 30;
        double A = l*sin(beta1);
        double B = l*cos(beta1);
        double C = (h + 0.5*D)*sin(beta1) - 0.5*D*tan(beta1);
        double E = (h + 0.5*D)*cos(beta1) - 0.5*D;
        y = A*sin(_x)*cos(_x) + B*pow(sin(_x),2) - C*cos(_x) - E*sin(_x);
        return y;
    }
    double diff(double _x)
    {
        double y;
        double D = 30;
        double A = l*sin(beta1);
        double B = l*cos(beta1);
        double C = (h + 0.5*D)*sin(beta1) - 0.5*D*tan(beta1);
        double E = (h + 0.5*D)*cos(beta1) - 0.5*D;
        y = -A*pow(sin(_x),2) + A*pow(cos(_x),2) + 2*B*sin(_x)*cos(_x) + C*sin(_x) - E*cos(_x);
        return y;
    }
};

int main()
{
    Fun1 f1;
    Newton n1(33*pi/180, eps, 20, f1);
    double ans1 = n1.solve();
    std::cout << "root using Newton method in (a) is alpha = " << ans1*180/pi << "°, f(alpha) = " << (f1)(ans1) << std::endl;

    Fun2 f2;
    Newton n2(33*pi/180, eps, 20, f2);
    double ans2 = n2.solve();
    std::cout << "root using Newton method in (b) is alpha = " << ans2*180/pi << "°, f(alpha) = " << (f2)(ans2) << std::endl;

    Secant s1(32*pi/180, 33*pi/180, eps, 20, f1);
    double ans3 = s1.solve();
    std::cout << "root using Secant method in (c) is alpha = " << ans3*180/pi << "°, f(alpha) = " << (f1)(ans3) << std::endl;   

    Secant s2(150*pi/180, 170*pi/180, eps, 20, f1);
    double ans4 = s2.solve();
    std::cout << "root using Secant method with different initial value in (c) is alpha = " << ans4*180/pi << "°, f(alpha) = " << (f1)(ans4) << std::endl;   

    return 0;
}
